Potential Carbon Storage Facilities Near Import/Export Ports
Contents
Potential Carbon Storage Facilities Near Import/Export Ports¶
Import and Procedural Functions¶
import pandas as pd
import matplotlib.pyplot as plt
import geopandas as gpd
import folium
import contextily as cx
import rtree
from zlib import crc32
import hashlib
from shapely.geometry import Point, LineString, Polygon
import numpy as np
from scipy.spatial import cKDTree
from shapely.geometry import Point
from haversine import Unit
from geopy.distance import distance
/Users/jnapolitano/venvs/finance/lib/python3.9/site-packages/geopandas/_compat.py:111: UserWarning: The Shapely GEOS version (3.10.2-CAPI-1.16.0) is incompatible with the GEOS version PyGEOS was compiled with (3.10.1-CAPI-1.16.0). Conversions between both will be slow.
warnings.warn(
Restrictions¶
Must be near a pipeline terminal
Must be Near a petrolium Port
Query Plan¶
Imports
Import LNG terminal Dataa
Import well data
Filtering
for each well calculate nearest terminal
for each well calculate distance from nearest terminal
eliminate wells further than 15 miles from a terminal
Data Frame Import¶
Wells DataFrame¶
## Importing our DataFrames
gisfilepath = "/Users/jnapolitano/Projects/data/energy/non-active-wells.geojson"
wells_df = gpd.read_file(gisfilepath)
wells_df = wells_df.to_crs(epsg=3857)
Terminal DataFrame¶
## Importing our DataFrames
gisfilepath = "/Users/jnapolitano/Projects/data/energy/Liquified_Natural_Gas_Import_Exports_and_Terminals.geojson"
terminal_df = gpd.read_file(gisfilepath)
terminal_df = terminal_df.to_crs(epsg=3857)
Eliminating SUSPENDED Terminal from the DataFrame¶
terminal_df.drop(terminal_df[terminal_df['STATUS'] == 'SUSPENDED'].index, inplace = True)
terminal_df.rename(columns={"NAME": "TERMINAL_NAME"})
terminal_df['TERMINAL_GEO'] = terminal_df['geometry'].copy()
terminal_df.columns
Index(['OBJECTID', 'TERMID', 'NAME', 'ADDRESS', 'CITY', 'STATE', 'ZIP', 'ZIP4',
'TELEPHONE', 'TYPE', 'STATUS', 'POPULATION', 'COUNTY', 'COUNTYFIPS',
'COUNTRY', 'LATITUDE', 'LONGITUDE', 'NAICS_CODE', 'NAICS_DESC',
'SOURCE', 'SOURCEDATE', 'VAL_METHOD', 'VAL_DATE', 'WEBSITE', 'EPA_ID',
'ALT_NAME', 'OWNER', 'POSREL', 'JRSDCTN', 'CONTYPE', 'IE_PORT',
'BERTHS', 'STORAGE', 'STORCAP', 'CURRENTCAP', 'APPCAP', 'OPYEAR',
'IMPEXPCTRY', 'VOLUME', 'PRICE', 'geometry', 'TERMINAL_GEO'],
dtype='object')
Plotting Terminal by TYPE¶
terminal_map =terminal_df.explore(
column="TYPE", # make choropleth based on "PORT_NAME" column
popup=True, # show all values in popup (on click)
tiles="Stamen Terrain", # use "CartoDB positron" tiles
cmap='Set1', # use "Set1" matplotlib colormap
#style_kwds=dict(color="black"),
marker_kwds= dict(radius=6),
#tooltip=['NAICS_DESC','REGION', 'COMMODITY' ],
legend =True, # use black outline)
categorical=True,
)
terminal_map
Terminal Impressions¶
According to the data there is not an export nor import location on The Western Side of the United States.
East Asian import or carbon capture export demands could justfity port development. Another study must be conducted.
Filtering Wells by Terminal Distance in Scipy¶
Edit¶
This method does not accuraletly calculate distance. The algorith used below calculates distance on a euclidan plane. In order to calculate a correct answer we must account for sphericity.
I include the code below as reference and a learning opportunity
def ckdnearest(gdA, gdB):
nA = np.array(list(gdA.geometry.apply(lambda x: (x.x, x.y))))
nB = np.array(list(gdB.geometry.apply(lambda x: (x.x, x.y))))
btree = cKDTree(nB)
dist, idx = btree.query(nA, k=1)
gdB_nearest = gdB.iloc[idx].drop(columns="geometry").reset_index(drop=True)
gdf = pd.concat(
[
gdA.reset_index(drop=True),
gdB_nearest,
pd.Series(dist, name='dist')
],
axis=1)
return gdf
c = ckdnearest(wells_df, terminal_df)
c.describe()
| index | OBJECTID | LATITUDE | LONGITUDE | PERMITNO | OPERATORID | SURF_LAT | SURF_LONG | BOT_LAT | BOT_LONG | ... | LATITUDE | LONGITUDE | BERTHS | STORAGE | STORCAP | CURRENTCAP | APPCAP | VOLUME | PRICE | dist | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 160432.000000 | 160446.000000 | ... | 160446.000000 | 160446.000000 | 160446.00000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 1.604460e+05 |
| mean | 3.596509e+05 | 3.596519e+05 | 37.195605 | -89.030287 | 1.980508e+06 | 6.624020e+07 | -63.243732 | -176.705395 | -949.703169 | -955.920831 | ... | 34.083355 | -87.268514 | -33.95252 | -30.998635 | -125.730124 | 1.641901 | -918.544126 | -131.529641 | -131.638380 | 7.437940e+05 |
| std | 2.806046e+05 | 2.806046e+05 | 4.464367 | 8.648121 | 3.210840e+06 | 2.528493e+08 | 306.472981 | 269.412702 | 220.324332 | 192.546614 | ... | 4.011683 | 12.477050 | 185.97266 | 186.552966 | 346.868737 | 0.518523 | 272.088477 | 344.706936 | 344.554868 | 6.029022e+05 |
| min | 0.000000e+00 | 1.000000e+00 | 26.047093 | -123.342380 | -9.990000e+02 | -9.990000e+02 | -999.000000 | -999.000000 | -999.000000 | -999.000000 | ... | 27.889869 | -116.847415 | -999.00000 | -999.000000 | -999.000000 | 0.000000 | -999.000000 | -999.000000 | -999.000000 | 5.272348e+02 |
| 25% | 1.392022e+05 | 1.392032e+05 | 32.810740 | -93.942380 | -9.990000e+02 | -9.990000e+02 | 31.861820 | -94.036088 | -999.000000 | -999.000000 | ... | 30.112960 | -93.288224 | 2.00000 | 2.000000 | 7.500000 | 1.500000 | -999.000000 | 0.000000 | 0.000000 | 3.672892e+05 |
| 50% | 3.128155e+05 | 3.128165e+05 | 38.293375 | -87.987450 | -9.990000e+02 | -9.990000e+02 | 38.067251 | -89.171775 | -999.000000 | -999.000000 | ... | 31.988522 | -88.503111 | 2.00000 | 4.000000 | 10.400000 | 1.800000 | -999.000000 | 0.000000 | 0.000000 | 5.391763e+05 |
| 75% | 6.171978e+05 | 6.171988e+05 | 39.201190 | -81.236495 | 3.305543e+06 | -9.990000e+02 | 39.186490 | -81.239420 | -999.000000 | -999.000000 | ... | 38.389603 | -76.409092 | 2.00000 | 7.000000 | 14.800000 | 1.800000 | -999.000000 | 9.100000 | 9.370000 | 9.642862e+05 |
| max | 1.506196e+06 | 1.506197e+06 | 50.744220 | -75.593800 | 2.017004e+07 | 1.044775e+09 | 50.744220 | -75.593800 | 45.179110 | -78.423180 | ... | 42.392856 | -71.058151 | 2.00000 | 7.000000 | 18.000000 | 4.000000 | 4.000000 | 932.000000 | 11.340000 | 3.191132e+06 |
8 rows × 24 columns
nearest_wells_df= wells_df.sjoin_nearest(terminal_df, distance_col="distance_euclidian")
nearest_wells_df.describe()
| index | OBJECTID_left | LATITUDE_left | LONGITUDE_left | PERMITNO | OPERATORID | SURF_LAT | SURF_LONG | BOT_LAT | BOT_LONG | ... | LATITUDE_right | LONGITUDE_right | BERTHS | STORAGE | STORCAP | CURRENTCAP | APPCAP | VOLUME | PRICE | distance_euclidian | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 160432.000000 | 160446.000000 | ... | 160446.000000 | 160446.000000 | 160446.00000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 1.604460e+05 |
| mean | 3.596509e+05 | 3.596519e+05 | 37.195605 | -89.030287 | 1.980508e+06 | 6.624020e+07 | -63.243732 | -176.705395 | -949.703169 | -955.920831 | ... | 34.083355 | -87.268514 | -33.95252 | -30.998635 | -125.730124 | 1.641901 | -918.544126 | -131.529641 | -131.638380 | 7.437940e+05 |
| std | 2.806046e+05 | 2.806046e+05 | 4.464367 | 8.648121 | 3.210840e+06 | 2.528493e+08 | 306.472981 | 269.412702 | 220.324332 | 192.546614 | ... | 4.011683 | 12.477050 | 185.97266 | 186.552966 | 346.868737 | 0.518523 | 272.088477 | 344.706936 | 344.554868 | 6.029022e+05 |
| min | 0.000000e+00 | 1.000000e+00 | 26.047093 | -123.342380 | -9.990000e+02 | -9.990000e+02 | -999.000000 | -999.000000 | -999.000000 | -999.000000 | ... | 27.889869 | -116.847415 | -999.00000 | -999.000000 | -999.000000 | 0.000000 | -999.000000 | -999.000000 | -999.000000 | 5.272348e+02 |
| 25% | 1.392022e+05 | 1.392032e+05 | 32.810740 | -93.942380 | -9.990000e+02 | -9.990000e+02 | 31.861820 | -94.036088 | -999.000000 | -999.000000 | ... | 30.112960 | -93.288224 | 2.00000 | 2.000000 | 7.500000 | 1.500000 | -999.000000 | 0.000000 | 0.000000 | 3.672892e+05 |
| 50% | 3.128155e+05 | 3.128165e+05 | 38.293375 | -87.987450 | -9.990000e+02 | -9.990000e+02 | 38.067251 | -89.171775 | -999.000000 | -999.000000 | ... | 31.988522 | -88.503111 | 2.00000 | 4.000000 | 10.400000 | 1.800000 | -999.000000 | 0.000000 | 0.000000 | 5.391763e+05 |
| 75% | 6.171978e+05 | 6.171988e+05 | 39.201190 | -81.236495 | 3.305543e+06 | -9.990000e+02 | 39.186490 | -81.239420 | -999.000000 | -999.000000 | ... | 38.389603 | -76.409092 | 2.00000 | 7.000000 | 14.800000 | 1.800000 | -999.000000 | 9.100000 | 9.370000 | 9.642862e+05 |
| max | 1.506196e+06 | 1.506197e+06 | 50.744220 | -75.593800 | 2.017004e+07 | 1.044775e+09 | 50.744220 | -75.593800 | 45.179110 | -78.423180 | ... | 42.392856 | -71.058151 | 2.00000 | 7.000000 | 18.000000 | 4.000000 | 4.000000 | 932.000000 | 11.340000 | 3.191132e+06 |
8 rows × 25 columns
Calculating Distance in Kilometers from Import/Export Terminal¶
#df.geopy.distance.distance(coords_1, coords_2).km
#df.apply(lambda row: distance(row['point'], row['point_next']).km if row['point_next'] is not None else float('nan'), axis=1)
# Thanks to https://stackoverflow.com/questions/55909305/using-geopy-in-a-dataframe-to-get-distances
nearest_wells_df['true_distance_km'] = nearest_wells_df.apply(lambda row: distance((row['LATITUDE_left'], row['LONGITUDE_left']), (row['LATITUDE_right'], row['LONGITUDE_right'])).km if row['geometry'] is not None else float('nan'), axis=1)
nearest_wells_df.describe()
| index | OBJECTID_left | LATITUDE_left | LONGITUDE_left | PERMITNO | OPERATORID | SURF_LAT | SURF_LONG | BOT_LAT | BOT_LONG | ... | LONGITUDE_right | BERTHS | STORAGE | STORCAP | CURRENTCAP | APPCAP | VOLUME | PRICE | distance_euclidian | true_distance_km | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 1.604460e+05 | 1.604460e+05 | 160446.000000 | 160446.000000 | 160432.000000 | 160446.000000 | ... | 160446.000000 | 160446.00000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 160446.000000 | 1.604460e+05 | 160446.000000 |
| mean | 3.596509e+05 | 3.596519e+05 | 37.195605 | -89.030287 | 1.980508e+06 | 6.624020e+07 | -63.243732 | -176.705395 | -949.703169 | -955.920831 | ... | -87.268514 | -33.95252 | -30.998635 | -125.730124 | 1.641901 | -918.544126 | -131.529641 | -131.638380 | 7.437940e+05 | 592.941775 |
| std | 2.806046e+05 | 2.806046e+05 | 4.464367 | 8.648121 | 3.210840e+06 | 2.528493e+08 | 306.472981 | 269.412702 | 220.324332 | 192.546614 | ... | 12.477050 | 185.97266 | 186.552966 | 346.868737 | 0.518523 | 272.088477 | 344.706936 | 344.554868 | 6.029022e+05 | 467.343457 |
| min | 0.000000e+00 | 1.000000e+00 | 26.047093 | -123.342380 | -9.990000e+02 | -9.990000e+02 | -999.000000 | -999.000000 | -999.000000 | -999.000000 | ... | -116.847415 | -999.00000 | -999.000000 | -999.000000 | 0.000000 | -999.000000 | -999.000000 | -999.000000 | 5.272348e+02 | 0.454603 |
| 25% | 1.392022e+05 | 1.392032e+05 | 32.810740 | -93.942380 | -9.990000e+02 | -9.990000e+02 | 31.861820 | -94.036088 | -999.000000 | -999.000000 | ... | -93.288224 | 2.00000 | 2.000000 | 7.500000 | 1.500000 | -999.000000 | 0.000000 | 0.000000 | 3.672892e+05 | 311.703128 |
| 50% | 3.128155e+05 | 3.128165e+05 | 38.293375 | -87.987450 | -9.990000e+02 | -9.990000e+02 | 38.067251 | -89.171775 | -999.000000 | -999.000000 | ... | -88.503111 | 2.00000 | 4.000000 | 10.400000 | 1.800000 | -999.000000 | 0.000000 | 0.000000 | 5.391763e+05 | 421.824879 |
| 75% | 6.171978e+05 | 6.171988e+05 | 39.201190 | -81.236495 | 3.305543e+06 | -9.990000e+02 | 39.186490 | -81.239420 | -999.000000 | -999.000000 | ... | -76.409092 | 2.00000 | 7.000000 | 14.800000 | 1.800000 | -999.000000 | 9.100000 | 9.370000 | 9.642862e+05 | 799.595274 |
| max | 1.506196e+06 | 1.506197e+06 | 50.744220 | -75.593800 | 2.017004e+07 | 1.044775e+09 | 50.744220 | -75.593800 | 45.179110 | -78.423180 | ... | -71.058151 | 2.00000 | 7.000000 | 18.000000 | 4.000000 | 4.000000 | 932.000000 | 11.340000 | 3.191132e+06 | 2390.537825 |
8 rows × 26 columns
Filtering Wells within 50 KM of a Terminal¶
filtered_wells = nearest_wells_df.loc[nearest_wells_df['true_distance_km'] < 50].copy()
filtered_wells.describe()
| index | OBJECTID_left | LATITUDE_left | LONGITUDE_left | PERMITNO | OPERATORID | SURF_LAT | SURF_LONG | BOT_LAT | BOT_LONG | ... | LONGITUDE_right | BERTHS | STORAGE | STORCAP | CURRENTCAP | APPCAP | VOLUME | PRICE | distance_euclidian | true_distance_km | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 2.654000e+03 | 2.654000e+03 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.0 | 2654.000000 | 2654.000000 | 2652.000000 | 2654.000000 | ... | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 | 2654.000000 |
| mean | 1.554318e+05 | 1.554328e+05 | 30.005777 | -93.585830 | 16514.931047 | -999.0 | 29.230172 | -94.272273 | -908.678451 | -919.670737 | ... | -93.524363 | -375.983798 | -375.082894 | -370.543293 | 0.973576 | -350.589940 | -374.631349 | -376.703112 | 24095.468770 | 20.848131 |
| std | 3.226352e+05 | 3.226352e+05 | 0.377108 | 0.665728 | 90735.288423 | 0.0 | 28.244840 | 24.858548 | 291.080631 | 255.761340 | ... | 0.630421 | 485.299749 | 486.001893 | 489.542071 | 1.090584 | 478.695054 | 488.693585 | 484.744059 | 14264.583731 | 12.326692 |
| min | 2.500000e+01 | 2.600000e+01 | 27.616626 | -97.747110 | -999.000000 | -999.0 | -999.000000 | -999.000000 | -999.000000 | -999.000000 | ... | -97.274403 | -999.000000 | -999.000000 | -999.000000 | 0.000000 | -999.000000 | -999.000000 | -999.000000 | 527.234827 | 0.454603 |
| 25% | 1.038925e+04 | 1.039025e+04 | 29.986090 | -93.621635 | -999.000000 | -999.0 | 29.985485 | -93.621915 | -999.000000 | -999.000000 | ... | -93.337457 | -999.000000 | -999.000000 | -999.000000 | 0.000000 | -999.000000 | -999.000000 | -999.000000 | 12302.653045 | 10.625148 |
| 50% | 5.860500e+04 | 5.860600e+04 | 30.032030 | -93.441460 | -999.000000 | -999.0 | 30.031910 | -93.442460 | -999.000000 | -999.000000 | ... | -93.336491 | 2.000000 | 3.000000 | 9.000000 | 0.710000 | 1.410000 | 0.000000 | 0.000000 | 26936.447255 | 23.256971 |
| 75% | 1.096650e+05 | 1.096660e+05 | 30.251400 | -93.338690 | -999.000000 | -999.0 | 30.251395 | -93.338887 | -999.000000 | -999.000000 | ... | -93.329114 | 2.000000 | 3.000000 | 11.000000 | 1.800000 | 4.000000 | 0.000000 | 0.000000 | 32607.124822 | 28.267032 |
| max | 1.503142e+06 | 1.503143e+06 | 30.561100 | -88.052220 | 833497.000000 | -999.0 | 30.561100 | -92.782800 | 30.077654 | -93.845118 | ... | -88.503111 | 2.000000 | 5.000000 | 18.000000 | 4.000000 | 4.000000 | 932.000000 | 11.340000 | 58104.197789 | 49.997185 |
8 rows × 26 columns
Map of Wells within 50 km of an Import/Export Terminal by Type¶
filtered_wells.explore(
column="STATUS_left", # make choropleth based on "PORT_NAME" column
popup=True, # show all values in popup (on click)
tiles="Stamen Terrain", # use "CartoDB positron" tiles
cmap='Set1', # use "Set1" matplotlib colormap
#style_kwds=dict(color="black"),
marker_kwds= dict(radius=6),
#tooltip=['NAICS_DESC','REGION', 'COMMODITY' ],
legend =True, # use black outline)
categorical=True,)